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The equations describing the first two terms of an asymptotic expansion of the 
solution of the planar turbulent mixing layer for values of the velocity ratio close 
to one are obtained. The first term of this expansion is the solution of the well- 
known time-evolving problem and the second, which includes the effects of the 
increase of the turbulence scales in the stream-wise direction, obeys a linear system 
of equations. Numerical solutions of these equations for a two-dimensional reacting 
mixing layer show that the correction to the time-evolving solution may explain the 
asymmetry of the entrainment and the differences in product generation observed 
in flip experiments. 


1. Introduction 

Time evolving simulations of the mixing layer are believed to capture many im- 
portant features of the dynamics of this flow, and are widely used because they 
are simpler to implement, less expensive, and less subject to uncertainties coming 
from approximate boundary conditions than the alternative space-evolving simula- 
tions. Some features of the real flow, however, are outside the framework of the 
time-evolving simulations. In particular, it is known that an incompressible mixing 
layer does not entrain equal amounts of fluid from each of the two free streams (Di- 
motakis, 1986 and references therein) and that, when the free streams carry passive 
scalars, the average composition in the molecularly mixed fluid in the mixing layer 
is nearer to the composition of the fast stream than to the composition of the slow 
stream. While it is not clear to what extent the second of these features is a conse- 
quence of the first or of some asymmetry of the process leading to molecular mixing 
inside the layer, none of the two can be captured by a time-evolving simulation, 
owing to the intrinsic symmetry of this artificial flow. 

The amount of molecular mixing and the variation of the mixed fluid mean con- 
centration across the layer depend on the Reynolds number, the degree of develop- 
ment of the layer, and the Prandtl or Schmidt number of the scalar; these factors 
determining in particular the extent to which free stream fluid is transported across 
the layer by the large scale motions before being molecularly mixed (Konrad 1976; 
Batt 1977; Mungal & Dimotakis 1984; Koochesfahani & Dimotakis 1986; Ivarasso 
& Mungal 1996), but the asymmetries mentioned above occur in any case. 
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The effect of the scalar field asymmetry is specially obvious when the species 
carried by the two streams are reactive and lead to a diffusion flame inside the mixing 
layer, in which a unit of mass of a stream (the fuel stream say) reacts with 5 units 
of mass of the other (the oxidizer stream). Then, flip experiments, consisting in 
exchanging the reactive species between the two streams and keeping constant their 
concentrations and all the other operating conditions (Mungal & Dimotakis 1984; 
Koochesfahani & Dimotakis 1986; Karasso & Mungal 1996), clearly show that the 
amount of product generated by the chemical reaction is larger when the reactant 
that is more consumed (the oxidizer if S > 1) is carried by the fast stream. 

To study the effect of these asymmetries, a correction to the time-evolving mixing 
layer formulation is worked out here using an asymptotic expansion for small differ- 
ences of the two stream velocities. The analysis follows the lines of Spalart (1986, 
1988). The new formulation is applied to a constant density two-dimensional reac- 
tive mixing layer with an infinitely fast, diffusion controlled reaction. In this limit, 
and assuming in addition that the Lewis numbers of the two reactants are equal to 
one, the chemical reaction can be easily accounted for following the evolution of a 
linear combination of the mass fractions of the reactants (the mixture fraction Z) 
which takes the values 0 and 1 in the oxidizer and fuel streams and is transported as 
a passive scalar. The mass fraction of the product (also with Lewis number equal to 
one) is a piecewise linear function of Z given by (6) below (see, e.g., Williams 1985). 

2. Formulation 

The turbulent mixing layer is a statistically stationary flow and, asymptotically, 
after an adjustment region, it is also statistically self-similar. Our aim here is to try 
to take advantage of these two properties to simplify the numerical computations. 

The large length and time scales of the turbulence are proportional to the stream- 
wise distance x in the self-similar state, and the variation of these scales is inex- 
tricably linked to the dynamics of the layer, being brought about by the process 
responsible for its evolution (vortex pairings in the classical view). This feature 
makes the numerical simulation costly because a long and wide stretch of the layer 
must be simulated in order to see the approach to its self-similar state. 

The opening angle of the mixing layer is often small. Formally this is true when 
the velocities of the two streams are close to each other, but in fact the angle is 
fairly small in any case. On one hand, this feature makes the numerical simulations 
even more expensive, because a very long stretch is needed to see the initial size of 
a typical vortex grow by a given factor. On the other hand, this feature is the basis 
of a well known simplification of the numerical treatment, based on the fact that 
the changes of time and length scales are small, and can be taken into account as a 
perturbation, over distances of the order of the size of one or a few vortices. 

This is so because for a layer growing by pairings of the large vortices, for which 
each vortex must undergo a number of these discrete events to approach the self- 
similar state, a small angle amounts to a small ratio of the size of a vortex to the 
distance it travels between successive pairings. But then, since the spacing of the 
large scale vortices is not much larger than their size, there are many neighboring 
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vortices to the left and to the right of a given vortex with sizes not very different 
from its own size, and it seems reasonable to think that the evolution of a vortex 
and of its nearest neighbors between successive pairings depends on the details of 
the evolution of only a number of other neighbors, whereas the effect of the rest of 
the mixing layer can be taken into account using a kind of mean field approximation 
without regard for individual features of the far vortices. 

Since the vortices move with a velocity intermediate between the velocities of 
the two streams, these conditions can be better put to use in a reference frame 
moving with that intermediate velocity. In this reference frame the evolution of 
a vortex can be described following the detailed dynamics of the flow over a span 
containing only a limited number of neighbor vortices. Here, using periodicity 
conditions in the stream- wise direction amounts, of course, to a temporal simulation 
of the flow, but the equations describing the small perturbations due to the growth 
of the turbulent scales can be easily worked out and solved along with those of 
the temporal simulation. It is worth recalling, however, that a definite error of 
a different type is always associated to the use of periodic conditions because the 
self-similarity and stationarity of the flow are statistical properties not applicable 
to individual realizations. Hopefully, the importance of this error decreases as the 
number of vortices included in the simulation increases. 

Restricting the computations to a finite span in the moving frame imposes a 
limitation on the time during which the evolution of the flow may be followed, 
owing to the growth of the vortices. In principle, a long span is required if the 
initial conditions are not close to the final self-similar state and a large number of 
pairings is necessary to approach that state. In practice, such a long span could 
perhaps be avoided carrying out computations on a shorter one for a moderate 
lapse of time and then replicating the results one or a few times in new adjacent 
spans, introducing appropriate phase shifts and a slight decrease of resolution before 
continuing the computation. 

Consider then a plane self-similar mixing layer between two incompressible streams 
with velocities U\ and U 2 < U\. Assume that e = (U\ — U 2 )/{U\ + U 2 ) 1 to en- 

sure that the angle of the layer, defined in any convenient way, is small of 0(e ) , 
though, as was mentioned before, this angle is probably small for any e. The aim 
is to follow the evolution of a few adjacent vortices, of initial characteristic size 6 
say, during a time of order t c = 26/ (U\ — U 2 ) corresponding to a few pairings. The 
growth of the mixing layer thickness is taken into account switching to the variables 
(x*,r),z*), where r 7 = y*/(ex*) and (x*,y*,z*) are the usual Cartesian coordinates. 
The statistical properties of the flow imply then that the x-averages coincide with 
the time averages and both are functions of rj only. Next, introducing a reference 
frame moving with velocity U m = (U\ + t/ 2 )/ 2 , the solution is sought in the form 

v* = U m (i + ev), p* = pU^P, with x* = U m t c (t + ex), ( 1 ) 

where the non-dimensional variables (v,p) are of order unity and the distances 
and time are scaled with 6 and t c , respectively. In terms of these variables the 
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Navier-Stokes equations and the equation for the mixture fraction are 

G,v, - - a] d v u = 0 
t -j- ex 

d t vi + = -G,p+ ^V 2 u, * , (2) 

d t Z + (v t G,)Z = -±—V 2 Z 

Kerr 

where G = (<9 r , 1 /(t + ex)d v , d z ) is a symbolic vector, v = (u,i; — r/( 1 + eu),u>), 
the Laplacian operator acting on each component of the velocity is V 2 = d xx + (1 + 
€ 2 7) 2 )/{t + ex) 2 + d zz - 2 ej]/(t + ex)d xl + 2 e 2 rj/{t + ex) 2 d n , Re = (U\ - U2)6/2v, 
and Pr is the Prandtl number. The solution of (2) can be sought as a power series 
in e, of the form (v, p, Z) = (v 0 , p 0 , Z 0 ) + e(vj, pj, Z x ) + • ♦ •. 

Carrying this expansion into (2) we find, at leading order, 

Go^o, — 0 

dtVo i + (vojGoj )t>o, = — Go.Po + 

dt Z q + {vqj Goj )Zq = IlePr^®^ 0 

where G 0 = (d x , l/td^, d z ), v 0 = (u 0 , vo - w o), and V§ = d xx + 1 /t 2 d^ + d zz . 

At the next higher order linear equations are obtained for (vi ,pi , Z \ ) having some 
forcing terms proportional to x (arising from the expansion of the denominators in 
(2)) and other that do not contain x explicitly. The solution of these equations 
is of the form (\ u p u Z i) = x (vi 0 ,Pio, Zio) + (v u , Pn , Z u ), where, as can be 
easily verified, (vjo, P io, Z\o) — 5<(v 0 ,po,2o) in order for the time averages to be 
independent of x, and {\\\,p\uZ n) satisfy 

7 ~j 

GoiVu- + dtUo — —dqUo = 0 
dtviii + (#o> Go } )vui + Gqj )v Qt — —Go t pn + — qi, + F t 
dtZ ii + (Oo, Go, )Z\\ + {vn } Go } )Zq + uod t Z q — = 

(v 2 0 Z n +2d It Z 0 - -jd xrt z^j 

where F = -uodtVo+T}Uo/tdqVo+2/Re(d It \o-r)/t d Tn v 0 )+(-d t p 0 +T)/td v po, 0, 0). 
The velocity at this order is therefore v = Vq + e(x 3 *v 0 -f Vu ), which can be written 
as Vo (x,rj,t + ex) + ev 11? and the pressure and mixture fraction are analogous. 
The first term admits a simple interpretation: as far as it is concerned the state of 
development of the flow is proportional to x, being slightly more evolved downstream 
of a point than it is upstream. Thus, even if periodic boundary conditions are used 
as an approximation for the leading order problem (3), as will be done in the 





2nd order corrections to time- evolving mixing layers 


211 


following section, this first term is not periodic. The effect of the second term will 
be discussed later. 

Applying x- and f-averages in suitable order to Eqs. (3) and to the continuity 


equation in (4) yields 


dv o 
dr} 

d , A du 0 


= 0 


d 

dr/ 


(vq+Po) = 0 


(5) 


d , — ,,-s dZ o 

— (ii 0 Zo) - rj—— = 0 
dr} dr} 

dv u du 0 


where the bars denote averaged variables. The first of these equations gives v 0 = 0 
(a constant v 0 has no effect on the dynamics; it can be set equal to zero by an 
0(e) change in the orientation of the x axis). The last equation in (5) implies 
that the mean normal velocity is of order e relative to the variation of the stream- 
wise velocity, as could have been expected for a region of O(e) aspect ratio. Using 
the second equation to eliminate yduo/dr/ and integrating, this equation yields 
Vn - = constant. Since = 0 for rj -» ±oo, the normal velocity i>n tends 

to the same constant value on both sides of the mixing layer, and this constant can 
be set equal to zero as for v 0 . Therefore the mixing layer does not introduce any 
perturbation in the free streams to this order and the ingestion of fluid by the layer 
is due only to the linear growth with x of its upper and lower apparent boundaries. 

Since the orientation of the ar-axis is well determined by the conditions v 0 — V\ = 0 
outside the mixing layer, the upper and lower boundaries can be defined on the 
basis of the usual thicknesses. For example, using the scaled momentum thickness 
A m = (Ui - u*)(u* - U 2 )dy* /[ex*(Ui -U 2 ) 2 ] and the scaled product thickness 

A p = Yp{Z) dy*/(ex*), where Y P (Z) is the product mass fraction given by the 

piecewise linear function 


f Z/Z a for 0 < Z < Z a 

p = \ (1 - Z)K 1 - Z a ) for Z a < Z < 1 


( 6 ) 


with Z a — 1/(1 + 5), the scaled upper and lower boundaries are 


i f±°° 1 f ±0 ° — e f 

Ai=± -*L = 1 (1 - T 2 i 

r i oo 

= ± Y p (Z)d V 

Jo 


« 0«11 dr] + 


( 7 ) 




and A m = A+ + A~, A p = A+ + A". 
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FIGURE 1. Sketch of the apparent boundaries and entrainment process. 


3, Results 

The two-dimensional forms of Eqs. (3) and (4) were numerically solved to find 
the effect of the first order corrections on the growth rate and the asymmetry of the 
layer. For this purpose the variable y = which is the non-dimensional normal 
distance divided by the factor 1 + ex/t, was used instead of y. The equations then 
take the form 

V • v 0 = 0 

d t \o + V • (v 0 v 0 ) = — Vp 0 + 4-V 2 v 0 

Re 

dfZo + V • (v 0 Zo) = ■■■■ ■ V 2 Zq 


( 8 ) 


V • v n = - d,u 0 } 
1 


dfv n + V • (v 0 v„ + v n v 0 ) = -Vpn + 4V 2 vn + <9,R 

Re 

dtZ\\ 4- V-(v 0 Zn + v„Z 0 ) = — - V 2 Zn 4- d t ( — uqZo 4- — — rr d x Zo 

RePr \ RePr 


) O) 


with R = — Mo v o + 2/Red x Vo — (po^ 0), and were solved with the boundary con- 
ditions v 0 T i = Vn — Z\\ — 0 and Zo = (1,0) for y — ► ±oo, and periodicity 
conditions in the stream-wise direction. The use of periodicity conditions is an 
approximation for which the only possible justification seems to be that in the 
present variables they are compatible with the spatial growth of the turbulence 
scales, and that, hopefully, they do not distort the solution too much if the period 
is sufficiently larger than the size of the vortices during most of the simulation. The 
initial conditions for the leading order variables were the hyperbolic tangent profiles 
i/o = tanh2y and Zq = |(1 4 tanh2y) plus perturbations proportional to the most 
unstable linear mode and one or two sub-harmonics with different amplitudes and 
phases. The variables Vn and Z\\ were initially zero. In the simulations Pr — 1 
and Re = 500-1000, based on the initial vorticity thickness. 

With these conditions, Vq, Zo — 1/2 and pn change sign under the transformation 
(x,y,t) — ► (-x,-y y t) while Vu, Z u , and p 0 are left invariant. The mixing layer 
grows symmetrically in first approximation, leading to A+ o = A“ 0 = A mo /2, 
whereas the correction to the growth rate is antisymmetric: A^ = — A” , — A, ni 
say (the same relations hold for the product thickness when 5 = 1). Here A mo and 
A mi are the slopes of straight lines fitted to 6 mo = ~ Jf^(l - ul)dy and S mi = 
— \ / 0 °° l dy. The fluxes crossing the upper and lower apparent boundaries of 
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FIGURE 2. Momentum thickness ( ) and upper and lower boundaries 

( an d ? respectively) for e = 0.45 (r = U 2 /U 1 = 0.38), Re = 500 

and Pr = 1. 

the layer, scaled with eU m x* , are (see Fig. 1) = (1 ± e)A± = [A mo ± e(A + 

2A mi ) 4- * • *]/2, where the momentum thickness is used for definiteness, and the 
entrainment ratio, defined here as E = <j>+ /<j>~, is 

E = l + e fl+2^-1 +0(e 2 ). (10) 

. ^mo . 

The total momentum thickness 6 m = <$ mo and the upper and lower apparent 
boundaries, <5+ = <5 mo /2 + e6 mi and 6~ = 6 m J 2 - e<5 m , respectively, are given in 
Fig. 2 for a representative case displaying a pairing. The numerical results show 
that the layer opens more towaid the slow stream ( A mi < 0), but this effect is over- 
balanced by the higher speed of the fast stream resulting in a E slightly greater 
than one. 

Figure 3 shows the product thickness for S = 1 and for S = 8 and 1/8. The 
last two values correspond to a flip experiment in which the fuel and the oxidizer 
streams are exchanged. The results show that the generation of product is higher 
when the reactant that is more consumed is carried by the fast stream. 

Both results are in qualitative agreement with the experimental data. Quantita- 
tive comparisons are meaningless given the two-dimensional character of the present 
simulations. 

The explanation of these results can be traced, of course, to the form of the 
forcing terms on the right hand sides of Eqs. (9). These terms depend only on 
time derivatives of the leading order solution, and the signs of some of them can 
be easily guessed. Thus, since the antisymmetric profile of uo(y) (the bar meaning 
here x-average) gets thicker with time, -d t u 0 in the continuity equation is more 



214 F. J. Higuera, J. Jimenez & A. Linan 

6 

5 
4 

5 „ * 

2 

1 

0 

0 10 20 30 40 50 

time 

Figure 3. Product thickness for 5 = 1 ( ; not affected by the asymmetry) 

and for 5 = 8 ( ) and 1/8 ( ), corresponding to a flip experiment, for 

e = 0.45, Re = 500 and Pr = 1 . : common value for 5 = 8 and 1/8 when 

6 = 0 . 

often positive than negative in the upper part of the layer, which amounts to a 
distribution of sources, and vice versa in the lower part, which amounts to sinks. 

For the same reason, — dfU^ should be positive, on average, everywhere in the 
layer, and the numerical results show that this is also true of the whole forcing term 
—dt(ul + p 0 ) + 2/Red xt uo in the ^-momentum equation. This amounts to a force 
pushing the fluid in the stream- wise direction and leading toaun predominantly 
positive. Therefore the average velocity Tio + tu\\ approaches its asymptotic value 
+ 1 in the upper stream faster than its asymptotic value —1 in the lower stream, 
which explains the results in Fig. 2. 

The forcing term in the equation for Z\\ has a complicated structure with bands 
of alternate signs . On average, however, it is positive (as could have been expected 
of the term — df(uoZo), due to the increase of thickness of the profile of Z o with 
time), leading Jo a Z\\ with a banded structure but predominantly positive. Hence 
the region of Z near 1 in the upper part of the layer is wider than the region of 
Z near 0 in the lower part. This provides an explanation for the results in Fig. 3 
because the upper region is responsible for a larger fraction of the product than 
the lower region when 5 is small and the average position of the flame is shifted 
toward the upper side, while the lower region is responsible for a larger fraction of 
the product than the upper region when 5 is large and the average position of the 
flame is shifted toward the lower side. 

Since the forcing terms depend only on the leading order solution, they could 
perhaps be evaluated from the results of a three-dimensional time-evolving simula- 
tion (as that of Roger &: Moser 1994), which would give indications on whether the 
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above observed trends hold also for that more realistic case. 

Finally it may be noted that the condition e <C 1, used here as the basis of a 
formal expansion, may not be necessary for some of the results to hold. As was 
mentioned before, the angle of the layer is moderately small for any value of e 
because the eddy turn over time is always shorter than the time between pairings, 
and this alone provides the required scale separation. 
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